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ABSTRACT 



,—l ' 13 February 2009 

pq: 

i-C ■ The triple asteroidal system (87) Sylvia is composed of a 280-km primary and two 

small moonlets named Romulus and Remus (Marchis et al 2005). Sylvia is located in 
q , the main asteroid belt, with semi-major axis of about 3.49 AU, eccentricity of 0.08 and 

?— i ■ 11° of orbital inclination The satellites are in nearly equatorial circular orbits around 

ryj \ the primary, with orbital radius of about 1,360 km (Romulus) and 710 km (Remus). 

. In the present work we study the stability of the satellites Romulus and Remus. In 

order to identify the effects and the contribution of each perturber we performed 
| numerical simulations considering a set of different systems. The results from the 3- 

i body problem, Sylvia-Romulus-Remus, show no significant variation of their orbital 

' elements. However, the inclinations of the satellites present a long period evolution 

7— I , with amplitude of about 20° when the Sun is included in the system. Such amplitude 

■ is amplified to more than 50° when Jupiter is included. These evolutions are very 

similar for both satellites. An analysis of these results show that Romulus and Remus 

. are librating in a secular resonance and their longitude of the nodes are locked to each 

■ other. Further simulations show that the amplitude of oscillation of the satellites' 
| inclination can reach higher values depending on the initial values of their longitude 
• of pericentre. In those cases the satellites get caught in an evection resonance with 

Jupiter, their eccentricities grow and they eventually collide with Sylvia. However, the 
. i-H [ orbital evolutions of the satellites became completely stable when the oblateness of 

■ Sylvia is included in the simulations. The value of Sylvia's J2 is about 0.17, which 

is very high. However, even just 0.1% of this value is enough to keep the satellite's 
orbital elements with no significant variation. 

Key words: asteroid, satellite, stability, oblateness 



1 INTRODUCTION 

The first triple asteroidal system discovered was (87) Sylvia 
(Marchis et al. 2005). The primary body is about 280 km, 
while the two small satellites, Romulus and Remus, are 
about 18 ± 4 km and 7 ± 2 km, respectively. For simplic- 
ity, in this paper we will call the primary of the system (87) 
Sylvia as "Sylvia". The satellites are in nearly equatorial cir- 
cular orbits around the primary, with orbital radius of about 
1,360 km (Romulus) and 710 km (Remus). In this paper we 
study the stability of the satellites Romulus and Remus. 
Sylvia is located in the main asteroid belt, with semi-major 
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axis of about 3.49 AU, eccentricity of 0.08 and 11° of or- 
bital inclination. Therefore, the two main perturbers of the 
system are the Sun and Jupiter. 

In order to identify the effects and the contributions 
of each relevant perturber of the system we explored several 
different dynamical systems. We studied the interactions be- 
tween the two satellites, the perturbation due to the Sun 
and due to the planet Jupiter. We identify two secular res- 
onances that play important roles on the dynamics of such 
systems. Finally, we investigate the contribution due to the 
oblateness of Sylvia. 

This paper has the following structure. In the next sec- 
tion we present our numerical simulations considering differ- 
ent dynamical systems, the 3-, 4- and 5-body problems. The 
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effects due to the interaction between Romulus and Remus 
are analysed in Section 3. The occurrence of the evection 
resonance between the satellites and Jupiter and its conse- 
quences are shown in section 4. Section 5 is devoted to the 
effect of the oblateness and finally, in Section 6, we present 
our conclusions. 



2 SIMULATIONS OF THE 3, 4 AND 5-BODY 
PROBLEMS 

We performed numerical simulations of several different dy- 
namical systems. In all cases we used the Gauss-radau inte- 
grator (Everhart 1985) and considered a time span of 5 x 10 4 
years, that corresponds to about 5 x 10 6 orbital periods of 
Romulus. Such time span was found to be suitable in order 
to capture all the main long term effects associated to the 
secular perturbations due to Jupiter. In order to check the 
accuracy of the numerical integrations we monitored the to- 
tal energy of the system. The relative error was always lower 
than 1CP 11 . The physical and orbital data of the system (87) 
Sylvia adopted in most of the simulations is given in Table 
1. 

First of all, we looked at the orbital evolution of Sylvia. 
Sylvia's orbit around the Sun, under Jupiter's perturbation, 
presents a very stable evolution (Figure 1). Its semi-major 
axis has only short period oscillations with small amplitudes, 
while the eccentricity and inclination show periodic secular 
evolution with amplitudes of Ae ~ 0.1 and A7 ~ 3°. That 
is the expected evolution due to the secular perturbation 
caused by Jupiter. Also in Figure 1 (in red) is shown the 
evolution of Sylvia's orbital elements under Jupiter's pertur- 
bation, considering only the secular terms of the three-body 
disturbing function (chapter 7 of Murray & Dermott 1999). 

Now we present the results of the numerical integrations 
for a set of different dynamical systems. Figure 2 shows the 
temporal evolution of the eccentricity (e), the inclination 
(7) and the longitude of the ascending node (fl) of Romu- 
lus (left column) and Remus (right column). In each figure 
the results from three different systems are plotted: the 3- 
body problem, Sylvia- Romulus-Remus (in red); the 4-body 
problem, Sylvia-Romulus- Remus-Sun (in green); the 5-body 
problem, Sylvia-Romulus- Remus-Sun- Jupiter (in blue). 

The results from the 3-body problem simulation, Sylvia- 
Romulus-Remus, show no significant variation of the satel- 
lite's orbital elements. The inclusions of the Sun and Jupiter 
do not affect the evolution of the eccentricities of the satel- 
lites. However, their orbital inclinations show a significant 
change. They present a periodic secular evolution, having 
amplitude of about 20° with period of almost 6 x 10 years 
(only Sun), and amplitude of about 58° with period of al- 
most 22 x 10 3 years (Sun and Jupiter). Their longitudes of 
the node either librate with amplitude of more than 120° 
(only Sun) or circulate (Sun + Jupiter). 

We also note that the evolution of i and f2 for both 
satellites present very similar amplitudes and frequencies. Of 
course, these evolutions are mainly due to the perturbations 
of the Sun (green) and Sun and Jupiter (blue), respectively. 
However, the fact of them being so similar is associated to 
a connection between Romulus and Remus, which will be 
discussed in the following section. 



3 THE ROMULUS-REMUS CONNECTION 

The gravitational interaction between Romulus and Remus 
produces different outcomes according to the dynamical sys- 
tem considered. In this section we concentrate the analysis 
on the temporal evolution of the satellites' orbital inclina- 
tions. We compare the results from our numerical integra- 
tions with results from the secular perturbation theory (see 
for instance Murray & Dermott (1999)). 

From the numerical integrations of the 3-body system, 
Sylvia-Romulus- Remus (Figure 2, in red), we see that Ro- 
mulus produces oscillations of small amplitude (< 4°) on Re- 
mus's inclination, which in return produces an even smaller 
amplitude of oscillation on Romulus's inclination (< 1°). 
These results are very close to those obtained from the sec- 
ular perturbation theory (Figure 3, first row). The direct 
coupling between Romulus and Remus, showing that when 
the inclination of one increases the other decreases, is clearly 
seen in the zoom of these plots (Figure 3, first row). 

The secular perturbation from the Sun on Romulus and 
on Remus, separately, produces oscillations with the same 
amplitude (~ 20°), but with different periods (Figure 3, 
second row). The period of oscillation for the inclination 
of Remus is more than twice that for Romulus. However, 
when the secular perturbation from the Sun and one of the 
satellites on the other satellite is computed, the amplitude 
and also the period of the oscillation of the inclination of 
both satellites are almost the same (Figure 3, third row). 
Actually, the gravitational interaction between Romulus and 
Remus is such that Remus's inclination follows very close the 
behaviour of Romulus's inclination. Such results are in very 
good agreement with our numerical integrations of the 4- 
body problem - Sylvia-Romulus-Remus-Sun (Figure 3, last 
row). 

This connection between Romulus and Remus is due 
to a secular resonance. The longitude of the ascending 
node of both satellites are librating (Figure 2, last row, 
in green). The evolution of the resonant angle, given by 
<j> = f^Rom — ^Rcm, shows a libration of very short ampli- 
tude (Figure 4, middle), i.e., the satellites' longitudes of the 
ascending nodes get locked and librate with almost the same 
frequency, f^Rom ~ f^Rem. This result is also verified by the 
results from the secular perturbation theory (Figure 4, left). 

Despite of the fact that in the 5-body problem, Sylvia- 
Romulus-Remus-Sun- Jupiter, the longitude of the ascending 
node of both satellites are circulating (Figure 2, last row, in 
blue), that same secular resonance occurs (Figure 4, right). 
Therefore, the evolutions of the orbital inclination and lon- 
gitude of the node of the two satellites are very much similar 
in the 5-body system due to this secular resonance. If the 
gravitational influence of one of the satellites is not taken 
into account their orbital evolution are very different from 
each other (Figure 5). The amplitude of Remus's inclina- 
tion is smaller than that of Romulus and their eccentricities 
grow erratically. So, we conclude that the direct gravita- 
tional perturbation of one satellite on the other does not 
produce any significant variation on their orbital elements. 
However, when these satellites are under strong perturba- 
tions from the Sun and Jupiter, the two satellites get locked 
in a secular resonance such that the orbital evolution of one 
is very similar to the other. 
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Table 1. Physical and orbital data of Sylvia's system (Marchis et al. 2005). 





Mass (Kg) 


a 


e 


in 


n(°) 


«(°) 


/(°) 


Orbital Period 


Sylvia 
Romulus 
Remus 


1.4780 x 10 19 
3.6625 x 10 15 
2.1540 x 10 14 


3.49 AU 
1356 Km 
706 Km 


0.08 
0.001 
0.016 


10.855 
1.7 
2.0 


266,195 
273 
314 


73,342 
101 
97 


8,51412 
81.88 
12.695 


6.52 years 
3.65 days 
1.38 days 




1 2 3 4 5 

time (10 4 years) 



Figure 1. : Temporal evolution of the semi-major axis, eccentricity, inclination, argument of pericentre and longitude of the ascending 
node of Sylvia. These are the results from the numerical integration of the 5-body problem - Sylvia-Romulus-Remus-Sun-Jupiter (in 
black) and from secular perturbation theory - Sylvia-Sun- Jupiter (in red). 
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Figure 2. : Temporal evolution of the eccentricity, the inclination and the longitude of the ascending node of Romulus (left column) 
and Remus (right column). In each figure are plotted the results from three different numerical integrations. 1) In red, the 3-body 
problem - Sylvia-Romulus-Remus. 2) In green, the 4-body problem - Sylvia-Romulus-Remus-Sun. 3) In blue, the 5-body problem - 
Sylvia- Romulus- Remus-Sun- Jupiter. 



4 THE EVECTION RESONANCE AND 
COLLISION 

The amplitude of oscillation of the satellites' inclinations 
can reach higher values depending on the initial value of 
Jupiter's longitude of pericentre. In Figure 6 we present 
the temporal evolution of the eccentricity and inclination 
of Romulus (left column) and Remus (right column), con- 
sidering different initial values of Jupiter's longitude of the 



pericentre. These are the results from the numerical inte- 
grations for the 5-body system Sylvia-Romulus- Remus-Sun- 



Jupiter, for ro.jup = 315° (black), c^jup 



^ju P = 10° (red) and zu. 



Jup 



15° (blue) 



= 320° (green), 
We verified that 



for —40° < zu] up < 15° the inclinations increase to higher 
values and the eccentricities suddenly grow until the satel- 
lites collide with Sylvia 

In those cases the satellites get caught in an evection 
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12345 12345 
time (1 4 years) time (1 4 years) 

Figure 3. : Temporal evolution of the orbital inclinations of Romulus (left column) and Remus (right column). In the first row are 
presented the results from the secular perturbation theory for the 3-body system - Sylvia-Romulus-Remus. It is also shown a zoom of 
the plots in the top right corner of the respective figure. In the second row are the results from the secular perturbation theory for the 
3-body systems: Sylvia-Romulus-Sun (left) and Sylvia-Remus-Sun (right) In the third row are presented the results from the secular 
perturbation theory for the 4-body system - Sylvia-Romulus-Remus-Sun. In the last row are presented the results from the numerical 
integration for the 4-body system - Sylvia-Romulus-Remus-Sun (reproduced from Figure 2, in green) 
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Figure 4. : Temporal evolution of the resonant angle, given by the difference between the longitudes of the ascending nodes of Romulus 
and Remus, <p = ^nom — ^Rem- The left plot is the result from the secular perturbation theory for the 4-body system - Sylvia- Romulus- 
Remus-Sun. The central plot is the result from the numerical integration for the 4-body problem - Sylvia-Romulus- Remus-Sun. The right 
plot is the result from the numerical integration for the 5-body problem - Sylvia-Romulus-Remus-Sun-Jupiter. 



resonance with Jupiter. The evection resonance usually oc- 
curs when the period of the longitude of the pericentre of 
the satellite is very close to the orbital period of the per- 
turber (Brouwer & Clemence, 1961, Vieira Neto et al, 2006, 
Yokoyama et al. 2008). In the present study the resonant 
angle is given by the difference between the satellite's lon- 
gitude of the pericentre and the mean longitude of Jupiter. 
In the first row of Figure 7 is presented an example of the 
evolutions of the evection angles of Romulus and Remus, for 
the case of initial zuj up = 90° . One can note that when the 
evection angle starts to librate the corresponding eccentric- 
ity starts to grow until reach a collision with Sylvia. 



5 OBLATENESS 

In this section the considered subject is the oblateness of 
Sylvia. First we discuss the limitations of the shape deter- 
mination of Sylvia from the observational images. Then, 
from secular perturbation theory (Ferraz-Mello, 1979), we 
analise the contribution of the Sun in comparison to that 
of Sylvia's oblateness on the satellite's orbital inclination. 
Finally, we present our numerical simulations for the 5- 
body system Sylvia-Romulus-Remus-Sun- Jupiter including 
the oblateness of Sylvia. 



5.1 Observational Constraints 

The J2 of Sylvia was first measured from the analysis of both 
satellite orbits by Marchis et al. (2005). The knowledge for 
each of them of its mean motion and semi-major axis con- 
strains simultaneously both the mass of Sylvia and its J2 
through the generalized Kepler's third law. The J2 measured 
was thereby derived to be 0.175 ± 0.050. The global shape 
of Sylvia was inferred from a convex inversion method of 
available ligthcurves by Kaasalainen et al. (2002). Further- 
more, Marchis et al. (2006) recorded with the Keck Adaptive 
Optics system a resolved image of Sylvia's primary at an an- 
gular resolution of 57 mas which agreed very well with this 
3D-shape model validating the convex model. Using the 3D- 
shape model and assuming an uniform mass distribution in 
the interior of the primary, we calculate a theoretical J2 of 
0.12, globally consistent with the measured value although 
slightly lower. The discrepancy may arise either from an in- 
terior heterogeneity or from the shape model itself which 



could have some local non-convexities, not taken into ac- 
count in the inversion process, able to further increase the 
J2 of Sylvia. 



5.2 Dynamical Effect 

In this section we consider the effect of Sylvia's oblateness 
(J2) on the dynamics of Sylvia's satellites. We define II = 
IexpiQ as a complex variable associated to the inclination 
(/) and longitude of the ascending node (SI). We use the 
indices 'o' for Romulus and 'e' for Remus and 'S' for the Sun. 
In what follows, we do not yet consider Sylvia's oblateness. 
Thus the secular equation for a satellite's II perturbed by 
the Sun is (Ferraz-Mello, 1979): 

dU a 
dt 
where 

3 Gins 



= -iA So Ilo + iA.s Ils, 



(1) 



4 n a'i 



(2) 



with an equivalent equation for Remus. Here, G is the grav- 
itational constant, m is a mass, n a mean motion, a is a 
semimajor axis and i is the imaginary unit. The solution of 
Eq. (1), in the approximation tig = constant, is: 



n o — K expig t + lis, 



(3) 



-A So . 

~ 2 year~ 1 , which corre- 



where K is a complex constant and g 

We obtain A So ~ 0.1108 x 10" 
sponds to a period of 5673 years, confirming the results of 
Section 3. For Remus, we find a period of 15012 years by 
the same method. The forced component lis induced by the 
Sun is thus responsible for the amplitude of the oscillation 
of I around 20°, since |IIs| ~ 10° and the initial conditions 
for I e and I are near 0°. 

We now add the second satellite to the secular equa- 
tions. To do this, we consider a complex vectorial variable 
II = (n o ,n e ). The equation for the secular variation of II 
(Ferraz-Mello, 1979) is similar to Eq. (1): 



dUg 
dt 



A II + i B U s 



where A is a 2 x 2 matrix with elements: 



A(l,l) = -(A e + A So ), 



(4) 



(5) 
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Figure 5. : Temporal evolution of the eccentricity, inclination, longitude of the ascending node and argument of pericentre of Romulus 
(left column) and Remus (right column). These are the results from the numerical integrations for the 4-body systems: Sylvia-Romulus- 
Sun- Jupiter (left) and Sylvia-Remus-Sun-Jupiter (right) 
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Figure 6. : Temporal evolution of the eccentricity and inclination of Romulus (left column) and Remus (right column). These are the 
results from the numerical integrations for the 5-body system Sylvia-Romulus- Remus-Sun- Jupiter, with different initial values of Jupiter's 
longitude of the pericentre: roj up = 315° (black); ro.j up = 320° (green); roj up = 10° (red); roj up = 15° (blue). 



A{l,2)=A e , 
A(2,1)=A , 
A{2,2) = -(A a +A Se ,) 
and B is a 2 x 1 vector with components: 

B(l)=A So , 
B(2) = A Sc , 



where 



1 Gm e L 
4 n a ai 
_ 1 Gm L 
4 n e ai 

r - —b w 

U ~ „2 U 3/2> 



(6) 
(7) 
(8) 



(9) 
(10) 



(11) 
(12) 
(13) 



where b^J 2 is a Laplace coefficient (Ferraz-Mello, 1979). The 
solution of Eq. (4) is quite similar to the one-dimensional 

case: 



n = N.E + C.U 



s. 



(14) 



where N is a 2 x 2 complex matrix whose rows are eigenvec- 
tors of matrix A and the corresponding eigenvalues are gi 



and (?2 in E = (exp igi, exp 132). C is a 'unitary' vector (1, 1). 
This results from the fact that A(j, 1) + A(j, 2) + B(j) = 0, 
for j = 1,2. As far as the forced component is concerned, 
both Romulus and Remus share a common one, given by 
the inclination of the Sylvia-Sun orbital plane with respect 
to Sylvia's equator. Thus the 20° amplitude of oscillation of 
7 e and I are the same as in the case of just one satellite 
(Fig. 3). As to the frequencies of oscillation, they are the 
eigenfrequencies of matrix A and we find them to be: 



gi = -0.108 x 10" 2 year _1 , 
32 = -0.832 x W~ 1 year~ 1 . 



(15) 
(16) 



These correspond to periods of around 75 and 5820 years, 
a long period forced by the Sun and a short period due to 
the satellites themselves. This short period is similar to the 
secular one defined by just the two satellites around Sylvia 
and the long period (similar to the period of Romulus in the 
Sylvia-Romulus- Sun problem) is shared by both satellites as 
shown in Fig. (3) (third and fourth row). 

We now introduce the oblateness of Sylvia defined by 
its J2. To better understand its effect we first include it in 
the secular equations for just one satellite and the Sun. In 
this way, Eq. (1) becomes: 
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Figure 7. : Temporal evolution of the evection angle, eccentricity and inclination of Romulus (left column) and Remus (right column). 
These are the results from the numerical integrations for the 5-body system Sylvia-Romulus- Remus-Sun- Jupiter, with initial value of the 
satellite's longitude of the pericentre equal to 90°. 



^ = -i(A So + A Jo ) + iA So Ils 



where 

A Jo = 



3 GJ-zm a yb 
2 



n ai 



(17) 



(18) 



Here, m sy is Sylvia's mass and b is Sylvia's equatorial radius. 
Now the solution of (17) is: 



n o = K exp ig t + cUs, 



(19) 



where K is a constant and: 



go = ~(A So +A Jo 



A So + A Jo 



(20) 
(21) 



The value of A Jo , with J 2 = 0.17 is 1.7401yeor _1 . Since 
Aj a 3> As , first I will oscillate with a much higher fre- 
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quency and second, the forced component, now multiplied 
by a very small constant c, will almost vanish. 

The equations for both satellites, the Sun and includ- 
ing Sylvia's oblateness will be like Eq (4), but the diagonal 
elements of matrix A are now: 



A(l,l) = -(A e + A So +Aj ), 
A(2,2) = -(A + A Se +A Je ). 



(22) 
(23) 



This will have important consequences on the eigenval- 
ues (proper frequencies) and the forced components. They 
now become: 



gi = —1.745 year 1 , 
g-2 = —17.04 year^ 1 , 
C = (6.35 x 1(T 4 , 2.75 x 10" s ), 



(24) 
(25) 
(26) 



where C = A^B. Thus the consequence of Sylvia's oblate- 
ness is to critically stabilize the satellite's orbits, since it 
introduces a frequency that is much faster than those gen- 
erated by the other gravitational perturbations. 



5.2.1 The effect of Jupiter 

The direct effect of Jupiter on Sylvia's satellites is negligible 
as compared with the Sun's effect. However, by considering 
Jupiter, Sylvia's orbital plane is no longer fixed but now 
precess with a frequency given by: 



A* = 7 _ 



1 G mj Lj 



?y a Sy 



(27) 



where the subscript 'J' refers to Jupiter and 'Sy' to Sylvia 
and 



Lj = 



aj 



3/2 



(28) 



Since Jupiter has a nonzero inclination with respect to 
the reference plane (ecliptic) Sylvia's inclination-node com- 
plex variable IT will be the sum of a constant forced compo- 
nent and a proper component precessing with frequency /i. 
Thus: 



Hsy = IIj + K exp —ifit 



(29) 



where K is a constant to be determined as far as one knows 
Ilsy for a specific time. 

Now to solve equations (1), (4) and (17), we must re- 
place the constant lis by Us v given by equation (29). The 
solution for each equation and a specific satellite can be gen- 
erally given by: 



n o = Il 0p + ki Kexp-ifit + k 2 Uj, 



(30) 



where H 0p is the proper component which will not change 
from the case without Jupiter, and 



ki = 
k 2 = 



A So +Aj-n 

A So 
A 8o +Aj' 



K, 



(31) 
(32) 



for the case of just one satellite here represented by Romulus 
(index 'o'). When both satellites are included then ki and 



k 2 are two dimension vectors given by: 
fci = -(A + _1 B, 



ki 



-A~ X B, 



(33) 
(34) 



where A and B are matrices already defined for the case 
without Jupiter and 'D' is the identity matrix. 

This powerfull effect of the oblateness was also explored 
in our numerical simulations. In Figure 8 we present the 
results from the numerical integrations for the 5-body sys- 
tem Sylvia-Romulus-Remus-Sun- Jupiter, with different val- 
ues for J2 of Sylvia: i) In purple, J2 = 0. ii) In blue, 
■h = 10~ 3 J2Syivia- Hi) In green, J 2 = 10~ 2 J2s y ivia- iv) In 
red, J2 = J2Syivia- The temporal evolution of the inclination 
and the longitude of the ascending node of Romulus (left 
column) and Remus (right column) clearly show that the 
value of J2 has to be much smaller (less than one thousand 
times) than the currently estimated value for Sylvia, in order 
to display any effect from other bodies' perturbation. 



6 CONCLUSIONS 

We have explored the dynamics of the satellites of aster- 
oid (87) Sylvia. When the oblateness of Sylvia is not taken 
into account, the satellites present a very interesting dynam- 
ical evolution. The perturbations from the Sun and Jupiter 
introduce a huge increase on the satellites orbital inclina- 
tions. Depending on the initial longitude of the pericentre 
of Jupiter, the satellites can be captured in a kind of evec- 
tion resonance. The longitudes of pericentre of the satellites 
get locked to the mean longitude of Jupiter. Such resonance 
forces the satellites eccentricities to grow exponentially and 
they eventually collide with Sylvia. It is also noted that the 
evolutions of Romulus and Remus are very similar. This is 
due to a secular resonance between them, which is caused by 
the perturbations from Sun and Jupiter. Finally, the com- 
plete stability of the satellites is guaranteed by the oblate- 
ness of Sylvia. We show that even just one thousand of the 
current value of J2 would be enough to keep the satellites in 
stable orbits. 

There are other triple asteroid systems in the main 
belt with similar characteristics, (45) Eugenia (Marchis et al 
2007) and (216) Kleopatra (Marchis et al 2008), which will 
be studied adopting a similar approach. 
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